Ideal Green Roof Locations: Melbourne

Authored by: Ryan Waites, Hannah Smith

Duration: 120 mins

Level: Beginner

Pre-requisite Skills: Python

Scenario

1. As a city planner, I want to identify the locations which could most benefit from retrofitting with a green roof.

2. As a building manager or owner of a residence, I want to visualise the potential of retrofitting my building for energy efficiency and greening.

Learning objectives

At the end of this use case you will be able to:

• work with spatial databases using Geopandas

• visualise spatial data on an interactive map

Why the interest in green roofs?

What is a green roof?

Green roofs are rooftops systems with plants in a growing medium, which are usually irrigated in drier climates. They can be:

• Extensive: growing in a shallower medium, generally low growing plants like succulents or grasses
• Intensive: growing in a deeper medium, with taller vegetation including shrubs or trees
• Semi-intensive: partway between

In Melbourne's climate even extensive green roofs may need to be irrigated.

How is this relevant to Melbourne?

Roofs make up about 23% of the total space in the City of Melbourne [1], and over 90% of rooftops which are suitable for extensive roofs are also suitable for intensive. Therefore, there is a significant potential for green roofs to contribute to the urban forest target the council has set for 2040 [2].

Building turnover in the centre of Melbourne is slow, therefore the best approach is to identify retrofitting opportunities in the suitable existing buildings [2].

Urban areas with an elevated temperature relative to the rural surroundings is known as the urban heat island (UHI) effect [3]. It can refer to both the urban surface temperatures on roads, footpaths, and building envelopes as well as the ambient temperature. There is evidence of the UHI in Melbourne with an average annual intensity of approximately 1.5 degrees Celsius [4].

The UHI threatens the economy by reducing productivity, increases the need for cooling and thus impacts energy consumption, and is associated with poorer air quality and a number of health risks [5]. In Melbourne, simulations of thermal environments for residential buildings indicate that doubling the amount of vegetation is estimated to reduce heat-related deaths by 5-28%. [6].

What are the benefits of green roofs?

Green roofs and walls mitigate the UHI by providing shade which blocks solar radiation from reaching urban surfaces, and by the vegetation absorbing and dissipating the radiation. This provides energy savings to the building owner, as well as cooling and comfort to the building users. For private businesses, this can improve productivity due to the improved interior comfort, and due to the psychological benefits of viewing vegetation [7].

In Melbourne, extensive green roofs were shown to reduce building energy use by 28% in summer [7]. Urban vegetation also increases urban ecology and biodiversity and provides amenity to the people using the urban spaces [8]. Green roofs make buildings cooler in summer and warmer in winter, provide barriers against noise pollution by reflecting sound as well as produce oxygen whilst capturing CO2 and other pollutants [9]. Ideally, trees should be used in combination with grasses, shrubs and planters for the optimum cooling effect [10], and it is more effective to use multiple UHI mitigation strategies at the same time (cool materials, green roofs, green walls, and urban greenery) [4].

In addition, green roofs are an effective and sustainable tool to assist with the management of stormwater quality and quantity. The rainwater is initially stored in the soil and then the vegetation, often reducing peak runoff and volume in contrast to conventional roofs. A study conducted in the highly urbanized city of Seoul, South Korea, found the green roof to have an average runoff retention ranging from 10% to 60% [11] whilst a study simulating Melbourne's rainfall patterns, commonly frequent and small, found they had the potential to retain up to 90% each rainfall [12]. These effects, in combination with the delaying of runoff, can greatly reduce the chances of flash flooding in highly urbanized areas and the load on urban drainage [13].

Relevant datasets

The Rooftop Project

This dataset uses a spatial multi-criteria analysis to classify buildings based on the possibility of adapting their rooftops for either intensive or extensive greenroofs. The buildings are located in the City of Melbourne.

• Rooftop datasets

---

Building Energy Consumption

This dataset outlines a model of energy consumption in megawatt hours in the City of Melbourne based on building attributes (age, floor area, etc) and is presented at property level scale. This model was developed by the CSIRO based on a baseline from 2011 to compare a potential reduction in energy consumption available in retrofitting is done. This retrofitting could be done demand-side (improving energy efficiency) or supply-side (by generating electricity on-site). The supply side considered roof size for placement of solar panels, but only considered 10% of that space to be potentially available for solar panels. The "business as usual" projections are property-level, but the retrofit scenario only is at block level due to privacy reasons.

• Property level energy consumption - business as usual
• Block level energy consumption - model for 2026

---

Urban Heat and Heat Vulnerability Index (HVI)

These datasets come from the Cooling and Greening Melbourne Interactive map made available by the state government. The urban heat dataset shows how many degrees Celsius the average temperature is above the baseline, when measuring within urban parts of a boundary area and using a non-urban baseline. The HVI refers to how vulnerable specific populations are to extreme heat events as indicated by heat exposure, sensitivity to heat, and adaptive capability.

• Access the Cooling and Greening Melbourne Interactive Map to download the data

---

The rooftop and urban heat/HVI datasets are not accessible via the Socrata API, they need to be downloaded from the links above onto your local machine. Please note that the Urban Heat and HVI datasets are only downloadable over a mixed security connection which will be automatically blocked by any Chromium browser.

# For importing from the Socrata API
import os
from sodapy import Socrata

# Working with the data
import numpy as np
import pandas as pd
import geopandas as gpd
from shapely.geometry import Polygon
from sklearn.preprocessing import LabelEncoder
import json

# Visualisation
import matplotlib.pylab as plt
import folium
from folium import plugins
import seaborn as sns
import branca.colormap as cm

%matplotlib inline

Accessing and loading the datasets

Loading the energy consumption datasets

The building energy consumption datasets are available by SODA API, available in the 'Export' tab at the website linked above. Since for this analysis the 2026 business as usual data will act as a baseline, it will be referred to as the baseline data.

# Importing using Socrata Open Data API using an anonymous app token

apptoken = os.environ.get("SODAPY_APPTOKEN")
domain = "data.melbourne.vic.gov.au"
client = Socrata(domain, apptoken)

WARNING:root:Requests made without an app_token will be subject to strict throttling limits.

# Import the datasets from the API endpoint and store as GeoPandas dataframes

baseline_dataset_id = "c5s6-vs5p"

baseline_df = gpd.GeoDataFrame.from_dict(
    client.get(baseline_dataset_id,
        # using SQL-like query to select only the columns we need
        select="the_geom, property_id, total_2026, floor_area",
        # if you don't specify a limit the dataset will be throttled to 1000 entries
        limit=15000))


model_dataset_id = "qucn-di27"

model_df = gpd.GeoDataFrame.from_dict(
    # This dataset has fewer than 1000 entries, so we don't need a limit
    client.get(model_dataset_id))

When downloaded from the API endpoint, these two datasets have their shapefiles stored as nested lists within a dictionary. To manipulate these datasets in Geopandas, we will need to convert these entries into a polygon.

Geopandas geometry columns are typically named 'geometry', so we will also rename the column. This avoids having to specify the column name when passing the shapefiles to Geopandas methods.

# Methods for preparing the energy consumption datasets

def get_coords(row):
    """
    Converts the geometry column entries of the datafile row into shapely polygons.
    row: a single row entry for a dataframe
    """

    row["geometry"] = Polygon(row['geometry']['coordinates'][0][0])

    return row


def prepare_dataset(dataframe, baseline=False):
    """
    Converts the downloaded energy consumption dataset contained in 'dataframe' into a usable Geopandas format
    using the 4326 coordinate system.

    Since we need different columns from each dataset, the bool 'baseline' will ensure the right columns are
    taken from each dataset.

    dataframe: a Pandas dataframe
    baseline: a bool indicating if this is the baseline or the model data
    """

    # Rename the column to be recognisable to Geopandas
    new_dataframe = dataframe.rename(columns={"the_geom" : "geometry"})

    # Apply the get_coords method to each row
    new_dataframe.apply(get_coords, axis=1)

    # The dataset also needs to have object datatypes converted to numeric data
    if baseline is True:
        new_dataframe["total_2026"] = pd.to_numeric(new_dataframe['total_2026'], errors = 'coerce')
        new_dataframe["floor_area"] = pd.to_numeric(new_dataframe['floor_area'], errors = 'coerce')
    if baseline is False:
        new_dataframe["total"] = pd.to_numeric(new_dataframe['total'], errors = 'coerce')

    return gpd.GeoDataFrame(new_dataframe,geometry='geometry', crs=4326)

# Preparing the energy consumption datasets

model_gdf = prepare_dataset(model_df,)
baseline_gdf = prepare_dataset(baseline_df, baseline=True)

Loading the rooftop datasets

Although the rooftop dataset is available as a zip over the API, this will not unpack to a shapefile. The shapefile data is only available by navigating to the link in the Relevant Datasets section and downloading the zipped file from there directly.

# Extensive green roof dataset
rooftop_ext_df = gpd.read_file("./interactive_dependencies/greenroof_usecase/rooftop_project/mga55_gda94_green_roof_extensive.shp")

# Intensive green roof dataset
rooftop_int_df = gpd.read_file("./interactive_dependencies/greenroof_usecase/rooftop_project/mga55_gda94_green_roof_intensive.shp")

# Solar roof dataset
rooftop_solar = gpd.read_file("./interactive_dependencies/greenroof_usecase/rooftop_project/mga55_gda95_green_roof_solar.shp")

# Cool roof dataset
rooftop_cool = gpd.read_file("./interactive_dependencies/greenroof_usecase/rooftop_project/mga55_gda94_green_roof_cool.shp")

Loading from the Cooling and Greening map

The urban heat dataset and the HVI dataset are downloaded by searching on the interactive map, which only allows downloading 1000 results at once. Four files are required to capture the entire area of the City of Melbourne.

Since the result selection on the interactive map is done with a rectangular selection, some of the data in the urban heat and HVI shape files will be irrelevant, so a shapefile of the City of Melbourne local government area (LGA) will be imported to trim the data to the correct areas.

# Heat vulnerability index dataset
hvi_df = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Heat Vulnerability Index (2018)(SA1).shp")

# The urban heat dataset
df_part1 = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Urban Heat (2018)(MB).shp")
df_part2 = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Urban Heat (2018)(MB)2.shp")
df_part3 = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Urban Heat (2018)(MB)3.shp")
df_part4 = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Urban Heat (2018)(MB)4.shp")
# Collating the three files
urbanheat_df = pd.concat([df_part1, df_part2, df_part3, df_part4]).drop_duplicates()

# The City of Melbourne LGA file
LGA_shape = gpd.read_file("./interactive_dependencies/greenroof_usecase/urban_heat/Local Government Area Solid.shp")

# Converting the datasets to the same coordinate system
rooftop_ext_gdf = rooftop_ext_df.to_crs(epsg=4326)
rooftop_int_gdf = rooftop_int_df.to_crs(epsg=4326)
rooftop_solar_gdf = rooftop_solar.to_crs(epsg=4326)
rooftop_cool_gdf = rooftop_cool.to_crs(epsg=4326)
hvi_gdf = hvi_df.to_crs(epsg=4326)
urbanheat_gdf = urbanheat_df.to_crs(epsg=4326)
LGA_shape_gdf = LGA_shape.to_crs(epsg=4326)

Filtering the data

This use case seeks to find which buildings in Melbourne City could most benefit from installing either an intensive or extensive green roof.

Filtering for the most suitable roofs

First, let's narrow down the green rooftop datasets to filter out any properties that aren't suitable. The RATING column uses categorical data, which can be encoded as numerical data to make it easier to work with. This converts the system from "Excellent" to "Very Poor" into 0-4.

# Converting categorical ratings to numerical form
le = LabelEncoder()

def label_encoded(column):
    # Encodes the data for the given column with labels between 0 and n_classes-1
    le.fit(column)
    le.classes_ = np.flip(le.classes_)
    return le.transform(column)

# Converting the columns in the dataset
# We don't need to make a copy of the data since this process is reversible using the .inverse_transform method
rooftop_ext_gdf["RATING"] = label_encoded(rooftop_ext_gdf["RATING"])
rooftop_int_gdf["RATING"] = label_encoded(rooftop_int_gdf["RATING"])
rooftop_solar_gdf["RATING"] = label_encoded(rooftop_solar_gdf["RATING"])
rooftop_cool_gdf["RATING"] = label_encoded(rooftop_cool_gdf["RATING"])

We need to visualise the spread of data across the rating categories to know how to filter. Below are bar charts which function as histograms showing the counts of each rating (ranging from "Very Poor" up to "Excellent") for both the intensive and extensive datasets.

Predicatably, the overall suitability for the intensive dataset is lower as this type of green roof requires more things, like suitability for irrigation.

# Styling for the plot
sns.set_theme(style="whitegrid", font_scale=1.2, palette="YlGn")

# Having the two plots share a y-axis allows direct comparison of the numbers between datasets,
fig, ((ax1, ax2), (ax3, ax4)) = plt.subplots(nrows=2,
    ncols=2,
    sharey=True, # both datasets will be scaled the same if they share a y-axis
    figsize=(10,8))

fig.tight_layout(pad=2.0)

sns.countplot(data=rooftop_ext_gdf, x=rooftop_ext_gdf["RATING"], ax=ax1)
sns.countplot(data=rooftop_int_gdf, x=rooftop_int_gdf["RATING"], ax=ax2)
sns.countplot(data=rooftop_solar_gdf, x=rooftop_solar_gdf["RATING"], ax=ax3)
sns.countplot(data=rooftop_cool_gdf, x=rooftop_cool_gdf["RATING"], ax=ax4)

# Shared properties for both axes
plt.setp((ax1, ax2, ax3, ax4),
    xticks=[4,3,2,1,0],
    xticklabels=['Excellent', 'Good','Moderate','Poor','Very Poor'],
    xlabel=" ",
    ylabel=" ")

# Properties for the extensive dataset
plt.sca(ax1)
plt.title("Extensive green roof suitability")
plt.ylabel("Count by rating category")

# Properties for the intensive dataset
plt.sca(ax2)
plt.title("Intensive green roof suitability")

# Properties for the extensive dataset
plt.sca(ax3)
plt.title("Solar roof suitability")
plt.ylabel("Count by rating category")

# Properties for the intensive dataset
plt.sca(ax4)
plt.title("Cool roof suitability")

plt.show()

There are several thousand properties which have an "excellent" suitability in each dataset, which is enough for the first pass of filtering.

# "Excellent" was encoded as 4 by the label encoder
extensive_gdf_best = rooftop_ext_gdf[rooftop_ext_gdf["RATING"] == 4]
intensive_gdf_best = rooftop_int_gdf[rooftop_int_gdf["RATING"] == 4]
solar_gdf_best = rooftop_solar_gdf[rooftop_solar_gdf["RATING"] == 4]
cool_gdf_best = rooftop_cool_gdf[rooftop_cool_gdf["RATING"] == 4]

These datasets have some very small polygons which might clutter up the map, so the final filtering for this dataset is to remove the very small shapes.

# Removing polygons with a small area
minimum_area = 100

extensive_gdf_best = extensive_gdf_best[extensive_gdf_best['Shape_Area'] > minimum_area]
intensive_gdf_best = intensive_gdf_best[intensive_gdf_best['Shape_Area'] > minimum_area]
solar_gdf_best = solar_gdf_best[solar_gdf_best['Shape_Area'] > minimum_area]
cool_gdf_best = cool_gdf_best[cool_gdf_best['Shape_Area'] > minimum_area]

Filtering the energy consumption datasets

Next, let's identify which buildings are modelled to have the smallest improvement in their energy efficiency.

To do this, we need to deal with the fact that the baseline dataset is property level, whereas the model dataset is block level.

Doing a spatial join on the datasets means that properties that share a block (as determined by the geometry column) will also share a right index. Since the total column is unique in the model_gdf dataset, that can be used to lookup the relevant row.

There are also zero values in both baseline_gdf and model_gdf for that need to be filtered out.

# Filtering the properties with totals of 0 mwh/s of energy consumption
model_gdf_filtered = model_gdf[model_gdf["total"] != 0]
baseline_gdf_filtered = baseline_gdf[baseline_gdf["total_2026"] != 0]

# Spatial joining the datasets
energy_consumption_diff = gpd.sjoin(baseline_gdf_filtered, model_gdf_filtered, how="inner")

We can visualise this by choosing a block to lookup, here I used block 358.

# Method to visualise the properties located on a given block in the energy consumption join

def visualise_properties(block_id):
    """
    Produces a plot visualising the properties located on a given block.

    block_id: the right_index of the block
    """
    fig, ax = plt.subplots(1, figsize = (5, 5))

    # Display the properties (baseline data)
    properties = energy_consumption_diff[energy_consumption_diff['index_right'] == block_id]

    # Display the block (modelling data)
    block = model_gdf_filtered[model_gdf_filtered['total'] == energy_consumption_diff[energy_consumption_diff['index_right'] == block_id].iloc[0]['total']]

    block.plot(ax=ax, color="#78c679")
    properties.plot(ax=ax, color="#31a354", lw=0.8)
    plt.yticks([])
    plt.xticks([])
    plt.title("Properties on block " + str(block_id))

    plt.show()

visualise_properties(358)

To find the difference between the projected business as usual condition for 2026 versus the modelled retrofitting scenario, we can create a new column and populate it with the percentage difference in energy consumption (increase or decrease). A positive number will represent the modelling situation having a greater energy consumption in megawatts per hour, and a negative number indicates a lower energy consumption.

We can then sort the dataset by the difference column.

To avoid redoing the summing operation for all the properties on a block, first sort the dataframe by the right index to ensure all properties in a block are clustered together.

# Adding a % difference between the baseline and model for 2026

# The new column to be added
difference_col = []

# Sort by the right join index
energy_consumption_diff = energy_consumption_diff.sort_values('index_right')

# Which index_right currently being worked on
current_index = 0
# The value to be stored in the current_index
current_val = 0

for row in energy_consumption_diff['index_right'].iloc:

    # If we are still processing the same index (block)
    if (current_index == row):

        # Add the same data as previous row
        difference_col.append(current_val)

    else:

        # Sum the baseline energy consumption for each property for this entire index
        baseline_sum = energy_consumption_diff.loc[energy_consumption_diff['index_right'] == row, 'total_2026'].sum()

        # Select the relevant total from the model total
        model = energy_consumption_diff.loc[energy_consumption_diff['index_right'] == row, 'total'].iloc[0]

        # Add that to the new difference column
        current_val = round((baseline_sum - model) / baseline_sum * -100, 2)
        difference_col.append(current_val)

        # Increment the current index so those calculations aren't redone unnecessarily
        current_index = row

# Add the new column and sort in descending order
energy_consumption_diff['difference'] = difference_col
energy_consumption_diff = energy_consumption_diff.sort_values('difference', ascending=False)

Inspecting the first few values of this data shows an unexpectedly large increase for some properties, given the difference column is a percentage. Block 180 has an increase of nearly 40000%.

Manual inspection of the blocks with unusually large increases reveals many of them have a floor area of 0, or few properties per block, which could be why the modelling total energy consumption is dramatically higher.

Filtering can remove the extraneously large increases in energy consumption.

visualise_properties(180)
visualise_properties(100)
visualise_properties(168)

# Selecting the top n properties in terms of increased energy consumption
# 1200 is around 10% of the total number of properties
n_properties = 1200

energy_consumption_diff_filtered = energy_consumption_diff[energy_consumption_diff['difference'] < 100]
energy_consumption_diff_filtered = energy_consumption_diff_filtered.iloc[:n_properties]

# Boxplot and histogram highlighting the distribution of values.
fig, ax = plt.subplots(1, 2, figsize = (16, 8))

energy_consumption_diff_filtered.boxplot(
    ax = ax[0],
    column = 'difference',
    patch_artist = True,
    boxprops=dict(facecolor = '#78c679',
    color = 'black'))


ax[0].set_title("Energy Consumption Differential Percentage")
ax[0].set_ylabel('Percentage')

energy_consumption_diff_filtered.hist(
    ax = ax[1],
    column = 'difference',
    color = '#78c679',
    edgecolor = 'black')

ax[1].set_title("Energy Consumption Differential Percentage")
ax[1].set_ylabel('Count')
ax[1].set_xlabel('Percentage')
ax[1].set_axisbelow(True)

energy_consumption_diff_filtered

	geometry	property_id	total_2026	floor_area	index_right	total	difference
5480	POLYGON ((144.97198 -37.81213, 144.97203 -37.8...	105969	97.701546	0.00	230	19579.580842	83.73
5478	POLYGON ((144.97141 -37.81289, 144.97131 -37.8...	105967	69.962541	0.00	230	19579.580842	83.73
3326	POLYGON ((144.97128 -37.81275, 144.97108 -37.8...	103625	1290.260295	0.00	230	19579.580842	83.73
3325	POLYGON ((144.97127 -37.81264, 144.97104 -37.8...	103624	1965.853186	0.00	230	19579.580842	83.73
1041	POLYGON ((144.97102 -37.81231, 144.97087 -37.8...	101117	274.619060	0.00	230	19579.580842	83.73
...	...	...	...	...	...	...	...
2369	POLYGON ((144.96690 -37.80674, 144.96691 -37.8...	102605	24.197633	204.84	303	4068.614405	-21.89
2383	POLYGON ((144.96689 -37.80559, 144.96714 -37.8...	102621	10.867908	92.00	303	4068.614405	-21.89
2379	POLYGON ((144.96685 -37.80578, 144.96739 -37.8...	102617	261.938489	426.87	303	4068.614405	-21.89
2373	POLYGON ((144.96676 -37.80639, 144.96678 -37.8...	102611	111.041668	940.00	303	4068.614405	-21.89
10296	POLYGON ((144.96687 -37.80621, 144.96732 -37.8...	111491	11.694814	99.00	303	4068.614405	-21.89

1200 rows × 7 columns

Filter the urban heat and HVI to the City of Melbourne LGA

The final filtering step concerns the urban heat difference and heat vulnerability index (HVI) datasets. Because the data was was selected by a rectangular area selection tool and we are only interested in the data that lies in the boundaries of the LGA, we need to filter out the values that don't relate to our area of interest - the City of Melbourne.

An intersection overlay with a Shapefile of the LGA boundaries will allow us to fit the datasets to the area we are interested in.

# Creating the overlays

HVI_overlay = gpd.overlay(hvi_gdf, LGA_shape_gdf, how="intersection")
urbanheat_overlay = gpd.overlay(urbanheat_gdf, LGA_shape_gdf, how="intersection")

# Visualising the overlay
fig, (ax1, ax2) = plt.subplots(1, 2, figsize = (16, 8))
ax1 = urbanheat_overlay.plot(ax=ax1, color="#31a354", lw=0.5)
ax2 = HVI_overlay.plot(ax=ax2, color="#31a354", lw=0.5)
ax1.set_axis_off()
ax2.set_axis_off()

# After filtering, examine the basic distribution of the dataset.
display(urbanheat_overlay['UHI18_M'].describe())

# Boxplot highlighting the distribution of values.
fig, ax = plt.subplots(1, 2, figsize = (16, 8))

boxplot = urbanheat_overlay.boxplot(
    ax = ax[0],
    column = 'UHI18_M',
    patch_artist = True,
    boxprops=dict(facecolor = '#78c679',
    color = 'black'))

ax[0].set_title("Urban Heat Index")

urbanheat_overlay.hist(ax = ax[1], column = 'UHI18_M', color = '#78c679', edgecolor = 'black')

ax[1].set_title("Urban Heat Index")
ax[1].set_axisbelow(True)

count    1537.000000
mean        7.983035
std         1.553946
min        -1.902073
25%         7.239289
50%         8.191541
75%         8.930228
max        12.019763
Name: UHI18_M, dtype: float64

# Sort the dataset from highest to lowest and examine the negative outliers
urbanheat_overlay.sort_values(by=['UHI18_M'], ascending = False, inplace = True)
display(urbanheat_overlay["UHI18_M"].tail(6))

# Map the lowest values against the rest of the area of interest
fig, ax = plt.subplots(1, figsize=(10, 10))
urbanheat_overlay.tail(6).plot(color='red', ax=ax)
model_gdf.plot(color='#78c679', ax=ax)
ax.set_axis_off()

1250   -0.762219
961    -0.835233
1251   -0.961165
1265   -1.454524
959    -1.767771
1252   -1.902073
Name: UHI18_M, dtype: float64

After a brief look at the distribution of values in the UHI for the Melbourne area, it's clear that the Urban Heat Island effect has a significant impact on temperatures in the city with an average increase of 8 degrees seen across all readings. Among all 1500+ readings, only 6 of manage to reach negative values and as to be expected, all occur in areas almost, if not entirely, covered by the river.

Visualisation

Finally, let's create a visualisation with our filtered data that allows us to locate buildings suitable for a green roof.

Since our rooftop datasets both have entirely "Excellent" properties, these can be visualised with a single colour. The other datasets will vary in colour based on a column's values.

The user will be able to interact with the map with each dataset appearing as a layer.

# Transform the encoded labels back to string for visualisation.
extensive_gdf_best["RATING"] = le.inverse_transform(extensive_gdf_best["RATING"])
intensive_gdf_best["RATING"] = le.inverse_transform(intensive_gdf_best["RATING"])
solar_gdf_best["RATING"] = le.inverse_transform(solar_gdf_best["RATING"])
cool_gdf_best["RATING"] = le.inverse_transform(cool_gdf_best["RATING"])

# Maps the rooftops datasets

def mapping(dataFrame, columns:list, labels:list, description: str, foliumMap):
    """
    dataFrame: the dataframe to add
    columns: list of the columns to be included in the map, with the most important column at index 0, 'geometry' must be last.
    labels: list of labels/aliases to use in place of provided column names for the tooltip (don't include one for 'geometry').
    description: a short description visible in map layers,
    foliumMap: the map to add the feature to
    """

    # Transform spatial data into a geoJson.
    geoJ = dataFrame[columns].copy()
    geoJ = geoJ.to_json()
    
    # Add the dataset name/description to the json properties as 'dataset' for use in tooltip.
    geoJ = json.loads(geoJ)
    for entry in geoJ['features']:
        entry['properties']['dataset'] = description
        
    # Add 'dataset' to the columns & labels lists for tooltip use.
    columns.insert(0, 'dataset')
    labels.insert(0, 'Dataset')
    
    # Feature group subgroup creation and adding to provided map (allows show/hide all selection).
    featureGroup = plugins.FeatureGroupSubGroup(rooftops_fg, name=description, overlay=True).add_to(foliumMap)

    # Add the features to the group
    feature = folium.GeoJson(
                geoJ,
                style_function=lambda feature: {
                    'color': 'black',
                    'fillOpacity': 0.2,
                    'opacity': .3,
                    'weight': 1
                    },
                # Controls behaviour when the user hovers the mouse over the feature
                highlight_function=lambda feature: {
                    'fillOpacity': 0.4,
                    'opacity': .5
                },
                # Hover tool tip will include all information except the 'geometry' column
                tooltip=folium.GeoJsonTooltip(fields = columns[:-1],
                                             aliases = labels)
                ).add_to(featureGroup)

# Mapping the data where the colours will vary based on a value

def mapping_overlay(dataFrame, columns:list, labels:list, description: str, foliumMap, colorMap, overlay):
    """
    dataFrame: the dataframe to add
    columns: list of the columns to be included in the map, with the most important column at index 0, 'geometry' must be last.
    labels: list of labels/aliases to use in place of provided column names for the tooltip (don't include one for 'geometry').
    description: a short description visible in map layers,
    foliumMap: the map to add the feature to
    colorMap: the colormap for this dataframe
    overlay: True or False value, sets as an overlay (Radio Button).
    """

    # Transform geo data into geoJson.
    geoJ = dataFrame[columns].copy()
    geoJ = geoJ.to_json()

    # Add the dataset name to json properties for use in tooltip.
    geoJ = json.loads(geoJ)
    for entry in geoJ['features']:
        entry['properties']['dataset'] = description
    
    # Add 'dataset' to the columns & labels lists for tooltip use.
    columns.insert(0, 'dataset')
    labels.insert(0, 'Dataset')
    
    # Feature group creation and adds the dataset to the map
    featureGroup = folium.FeatureGroup(name=description, overlay=overlay).add_to(foliumMap)

    # Add the features to the group
    feature = folium.GeoJson(
                geoJ,
                style_function=lambda feature: {
                    'fillColor': colorMap(feature['properties'][columns[1]]),
                    'fillOpacity': 0.5,
                    # outline for the overlayed datasets
                    'weight': 0.5,
                    'opacity': 0.5,
                    'color': 'black'
                    },
                # Controls behaviour when the user hovers the mouse over the feature
                highlight_function=lambda feature: {
                    'fillOpacity': 0.7,
                    'opacity': 0.7
                },
                # Hover tool tip will include all information except the 'geometry' column
                tooltip=folium.GeoJsonTooltip(fields=columns[:-1],
                                             aliases=labels)
                ).add_to(featureGroup)

# Create step color map for Urban Heat Index
heat_cmap = cm.StepColormap(
    colors=['gray', 'yellow', 'orange', 'orangered'],
    vmin = 0,
    vmax = 13,
    index = [0, 3, 6, 9, 12.2],
    caption = 'Urban Heat Index')

# Create step color map for Heat Vulnerability Index
vuln_cmap = cm.StepColormap(
    colors=['LightBlue', 'Plum', 'Orchid', 'MediumOrchid', 'DarkOrchid'],
    vmin = 0,
    vmax = 5,
    index =[0, 1, 2, 3, 4, 5],
    caption = 'Heat Vulnerability Index')

# Create base map
m = folium.Map(
    location=[-37.81368709240999, 144.95738102347036],
    zoom_start=13,
    width=1000,
    height=600,
    min_zoom=10,
    tiles=None)

# Adding the base 'always on' tile map (control on/off disabled to allow UHI/HVI overlay selection with base map still).
folium.TileLayer('Cartodb Positron', min_zoom=10, overlay=False, control=False).add_to(m)

# Add same as selectable layer to turn off all other tile layers.
folium.TileLayer('Cartodb Positron', min_zoom=10, overlay=False, name="No Overlay",).add_to(m)


mapping_overlay(
    urbanheat_overlay,
    # Can pass more columns, as long as the first column passed is the value to show in the color mapping
    ["UHI18_M", "geometry"],
    # Ensure a label for each column to be shown in the hover tool tip, except for 'geometry' which wont be listed.
    ["UHI in °C"],
    "Urban heat difference",
    m,
    heat_cmap,
    False)

mapping_overlay(
    HVI_overlay,
    ["HVI_INDEX", "geometry"],
    ["HVI"],
    "Heat Vulnerability Index",
    m,
    vuln_cmap,
    False)

mapping_overlay(
    energy_consumption_diff_filtered,
    ["difference", "floor_area", "geometry"],
    ["Difference %", "Total Floor Area"],
    "Modelled versus business as usual energy consumption",
    m,
    cm.linear.YlGn_09.scale(-100, 100),
    True)

# Create feature group to encompass all rooftop data (Enables show/hide all selection. Brings rooftops back above overlays.).
rooftops_fg = folium.FeatureGroup(name='Show/Hide Rooftops')
rooftops_fg.add_to(m)

mapping(
    extensive_gdf_best,
    ["RATING", "geometry"],
    ["Rating"],
    "Extensive Rooftop: Excellent",
    m)

mapping(
    intensive_gdf_best,
    ["RATING", "geometry"],
    ["Rating"],
    "Intensive Rooftop: Excellent",
    m)

mapping(
    solar_gdf_best,
    ["RATING", "geometry"],
    ["Rating"],
    "Solar Rooftop: Excellent",
    m)

mapping(
    cool_gdf_best,
    ["RATING", "geometry"],
    ["Rating"],
    "Cool Rooftop: Excellent",
    m)

folium.LayerControl().add_to(m)   # Add the layer control to switch layers

<folium.map.LayerControl at 0x1e51d11d220>

Finally, we output the map. Using the layer control button the overlays can be shown or hidden and the user can choose which of the filtered rooftop datasets to show at any one time.

# Displays the map
m

Make this Notebook Trusted to load map: File -> Trust Notebook

References

[1] GHD (2015) Rooftop Adaptation Study: Green Roofs, Cool Roofs and Solar Panels Final Report, City of Melbourne, Melbourne.

[2] Jones, R (2018) Valuing Green Guide. Green Roofs, Walls And Façades Project Report, City of Melbourne, Melbourne, Australia.

[3]: Howard L (1818) The Climate of London: Deduced from Metereological Observations, Made at Different Places in the Neighbourhood of the Metropolis, Howard W. Phillips, London, UK.

[4] Yenneti K, Ding L, Prasad D, Ulpiani G, Paolini R, Haddad S and Santamouris M (2020) 'Urban Overheating and Cooling Potential in Australia: An Evidence-Based Review'. Climate, 8(11):126. https://doi.org/10.3390/cli8110126

[5] Price A, Jones EC and Jefferson F (2015) 'Vertical Greenery Systems as a Strategy in Urban Heat Island Mitigation'. Water Air Soil Pollut 226, 247. https://doi.org/10.1007/s11270-015-2464-9

[6] Chen D, Wang X, Thatcher M, Barnett G, Kachenko A and Prince R (2014) 'Urban vegetation for reducing heat related mortality', Environmental Pollution 192:275-284.

[7] Williams N, Rayner J, Lee K, Fletcher T, Chen D, Szota C and Farrell C (2016) 'Developing Australian green roofs: Overview of a 5-year research program'. Acta Horticulturae. 1108:345-352. 10.17660/ActaHortic.2016.1108.46.

[8]: Hayes AT, Jandaghian Z, Lacasse MA, Gaur A, Lu H, Laouadi A, Ge H and Wang L (2022) 'Nature-based solutions (NBSs) to mitigate urban heat island (UHI) effects in Canadian cities', Buildings 12(925). https://doi.org/10.3390/buildings12070925

[9] Guattari C, Evangelisti L, Asdrubali F and De Lieto Vollaro R (2020) 'Experimental evaluation and numerical simulation of the thermal performance of a green roof', Applied Sciences 10(1767), https://doi.org/10.3390/app10051767

[10] Fu J, Dupre K, Tavares S, King, D and Banhalmi-Zakar Z (2022) 'Optimized greenery configuration to mitigate urban heat: A decade systematic review', Frontiers of Architectural Research 11:466-491.

[11] Jusić S, Hadžić E and Milišić H (2019) 'Stormwater Management by Green Roof', Acta Scientific: Agriculture, 3(7):57-62, DOI: 10.31080/ASAG.2019.03.0516

[12] Zhang Z, Szota C, Fletcher TD, Williams NSG and Farrell C (2019) 'Green roof storage capacity can be more important than evapotranspiration for retention performance', Journal of Environmental Management, 232:404-412, DOI: 10.1016/j.jenvman.2018.11.070

[13] Shafique M, Kim R and Kyung-Ho K (2018) 'Green Roof for Stormwater Management in a Highly Urbanized Area: The Case of Seoul, Korea', Sustainability, 10(3):584, DOI: 10.3390/su10030584

Return to top